### REPLICATION FILE -- ANALYSIS
### Jonathan Homola
### "Are Parties Equally Responsive to Women and Men?"
### British Journal of Political Science

## clear environment, set seed, install/load packages
rm(list=ls())
set.seed(12435); options(stringsAsFactors=F)
# install.packages("stargazer"); install.packages("sandwich")
# install.packages("psych"); install.packages("foreign")
library("foreign"); library("stargazer"); library(sandwich); library(psych)

## set working directory
setwd(" ... ")

## read in the main dataset 
data <- read.dta("Homola_WomenMenResponsiveness_Data.dta")

#### Table 2: Descriptive statistics: dependent and independent variables
data <- subset(data, niche==0)
data <- subset(data, huber==1)
data <- subset(data, !is.na(pshift1))
data <- subset(data, !is.na(mvshift))
data$absp <- abs(data$pshift)
data$absmv <- abs(data$mvshift) 
data$absw <- abs(data$wshift)
data$absm <- abs(data$mshift)
describe(data[,c("pshift", "mvshift", "wshift", "mshift", "absp", "absmv",
                 "absw", "absm", "womenpar")])[,c(2,3,4,8,9)]


#### Table 3: Parties' responsiveness to different electorates
mod1 <- lm(pshift ~ mvshift + pshift1, data)
mod2 <- lm(pshift ~ wshift + pshift1, data)
mod3 <- lm(pshift ~ mshift + pshift1, data)
mod4 <- lm(pshift ~ wshift + mshift + pshift1, data)

se1 <- sqrt(diag(vcovHC(mod1, "HC1")))
se2 <- sqrt(diag(vcovHC(mod2, "HC1")))
se3 <- sqrt(diag(vcovHC(mod3, "HC1")))
se4 <- sqrt(diag(vcovHC(mod4, "HC1")))

stargazer(mod1, mod2, mod3, mod4,
          se = list(se1, se2, se3, se4),
          omit.stat=c("f", "ser", "bic", "ll", "adj.rsq"), 
          star.cutoffs = c(0.1, 0.05, 0.01),
          column.sep.width="1pt", digits=2,
          dep.var.labels.include = FALSE,
          dep.var.caption="Outcome variable: Change in Party Position")


#### Table 4: Party responsiveness at different levels of female seat share
mod1 <- lm(pshift ~ wshift + womenpar + wparlxwshift + pshift1, data)
mod2 <- lm(pshift ~ wshift + womenpar + wparlxwshift + pshift1 + 
             votediff1 + pshiftxvotediff1, data)
mod3 <- lm(pshift ~ mshift + womenpar + wparlxmshift + pshift1, data)
mod4 <- lm(pshift ~ mshift + womenpar + wparlxmshift + pshift1 + 
             votediff1 + pshiftxvotediff1, data)
se1 <- sqrt(diag(vcovHC(mod1, "HC1")))
se2 <- sqrt(diag(vcovHC(mod2, "HC1")))
se3 <- sqrt(diag(vcovHC(mod3, "HC1")))
se4 <- sqrt(diag(vcovHC(mod4, "HC1")))

stargazer(mod1, mod2, mod3, mod4,
          se = list(se1, se2, se3, se4),
          omit.stat=c("f", "ser", "bic", "ll", "adj.rsq"), 
          star.cutoffs = c(0.1, 0.05, 0.01), 
          column.sep.width="1pt", digits=2,
          dep.var.labels.include = FALSE,
          dep.var.caption="Outcome variable: Change in Party Position")